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Abstract 


A method for designing robust feedback 
controllers for multi loop systems Is presented. 
Robustness Is characterized In terms of the 
minimum singular value of the system return 
difference matrix at the plant Input. Analytical 
gradients of the singular values with respect to 
design variables In the controller are derived. A 
cumulative measure of the singular values and 
their gradients with respect to the design 
variables Is used with a numerical optimization 
technique to Increase the system's robustness. 

Both unconstrained and constrained optimization 
techniques are evaluated. Numerical results are 
presented for a two-input/two-output drone flight 
control system. 

Nomenclature 


0 sideslip angle (deg) 

61,62 elevon and rudder deflections 

(deg) 

Of, nth singular value 

o*,o_ maximum and minimum singular value 

2 >i»£j) global minimum and desired 

singular value 
nth loop phase In L' matrix 

roll angle and rate (deg/sec) 

yaw angle and rate (deg/sec) 
uj frequency (rad/sec) 

[ ]* complex conjugate transpose of [ ] 

(*1 represents time derivative of ( ) 

tr[ ] trace of a square matrix [ ] 


Introduction 


A,B,C,D controller matrices 


augmented system matrices 
dB decibel 

F,Gy,H plant matrices 


G 

g 

I 

[ I+KG] 

J 

j 

K 

L 

M 

Ns.Nc.Nq 


P 

r 

s 

u 

Un.Vn 

X 


plant transfer matrix 
cumulative constraint 
Identity matrix 
return difference matrix 
elective function 

controller transfer matrix 

nth loop gain In L matrix 

diagonal gain and phase change matrix 

order of controller 

order of plant, Input, and 

output 

element of controller quadruple 
reference Input 
Laplace variable 
plant Input vector 
left and right eigenvectors 
plant state vector 
controller state vector 
plant output vector 


*Ae r os pa ce Engineer, HuTt Idlscipllnary Analysis 
and Optimization Branch, Loads and Aeroelasticity 
Division. 

tGeorge Washington University, Joint Institute for 
Advancement of Flight Sciences. 


A well -designed feedback control system 
should provide stability robustness with respect 
to plant uncertainty. For single-input/single- 
output systems, the classical concepts of gain and 
phase margins are employed as measures of system 
robustness. In multiloop systems, these classical 
single-loop measures may not always provide a good 
measure of system robustness. Recently, matrix 
singular value properties of a multiloop system's 
return difference matrix have been proposed as a 
measure of system robustness (1-4). Several 
authors have even related the singular values of 
the return difference matrix to multlloop gain and 
phase margins (3-4). 

The majority of the effort to date has 
focused on singular values as analysis tools. 

Only a small amount of effort has been focused on 
the use of singular values for control law 
synthesis (5-7). Stein (5) discusses the 
frequency domain Interpretation of the linear 
quadratic Gaussian (LQG) based design in terms of 
singular values. He shows how the LQG methodology 
can be used to design feedback controllers which 
satisfy design requirements expressed as singular 
value conditions. Safonov and Chen (7) discuss a 
procedure for maximizing singular values for 
stability margin optimization. The purpose of 
this paper Is to Introduce a new design method 
which employs a numerical optimization technique 
to search for the controller design variables that 
Increase the minimum singular value of the system 
return difference matrix. The singular value 




gradients required in the optimization schemes are 
derived analytically. Numerical results are 
computed for a two-input /two-output system which 
represents an experimental drone aircraft with a 
lateral attitude control system (8). 


System Description 


Let the multi loop feedback system shown in 
figure 1 be described by a set of constant 
coefficient differential equations of the form 

Plant 


^(L-l-I) ->(l-l/kn)i^ + 2/kn(l.cos*n) (8) 
max 
n 

' • *N^ 

Equation (8) is plotted in figure 2 with k^ and 
AS parameters. This figure can be used to 
determine the gain margins for a particular phase 
margin for simultaneous changes of both gain and 
phase in all input channels (4). 


Singular Value Gradient Derivation 


X = Fx + GyU 


z = Hx 


Controller 


Xc “ AX(^ + Bz 


^ ' In order to perform the optimization, it is 

necessary to determine the gradients of the 
/ 2 \ singular value o_(I+KG) with respect to elements in 

' ' the controller quadruple matrices A,8,C, and D. 

Let the parameter p represent one of the elements 
of the controller matrices which are the design 
variables. It was shown in reference 4 that for a 
distinct singular value 0 ^ of a complex matrix 
(3) (I+KG), the gradient with respect to a real 

parameter p is given by 


u ■ Cxc + Dz 


(4) 


Equation (1) represents an Ng order plant having 
Nq output measurements, z, modeled by equation 
(2) and N^ control inputs, u. Equations (3) and 
(4) represent an Mth order feedback controller 
driven by the sensor output z. In terms of a 
transfer function matrix, the plant and the 
controller are 


3oy^(UKG) 

ap 


Re (u„^ 


3(UKG) 

"ap 


Vfll 


(9) 


where v and u are respectively right and left 
normalized eigenvectors of (I+KG). (For repeated 
eigenvalues see reference 9 for the corresponding 
"Gateaux differential" expressions.) 

It can be shown that 


z = [H(Is-F)-l6u]u = G(s)u 

(5) 

u = [C(Is-A)‘^B + D]z = -K(s)z 

(6) 


respectively. 

Assuming the closed-loop system to be stable, 
the robustness of the nominal system at the plant 
input can be examined by computing o (I+KG) as a 
function of frequency (s-Jw) and usTng the 
guaranteed stability criterion 


o(L-l-I) <o(I+KG) (7) 

at all frequencies (3). In this paper, the matrix 
L Is a diagonal gain and phase change matrix at 
the Input of the plant as shown In figure 1. 


L ° Dlag [kne^^n] 


The matrix L Is the Identity matrix for the 
nominal system and It can be shown that 


3on(I+KG) 


Re [fi(Is-7f)-l^(vnUn*) 


[-i:r(Is-7I)-l7]] 


where 



(Nc+M)x(No+M) 


(10) 


(No+M)x(Nc+M) 



(No+M)x(Ns+M) 



(Ns+M)x(Ns+M) (Ns+M)xNc 


T - [ -0H| -Cl 

I 



Ncx(Ns+M) (Ns+M)xM 


N 
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To derive the matrix equation (10), define p 
as an element of the matrix P. Then the scalar 
equation (9) can be written as 


3on(UKG) 

TP 


Re* tr 


dCI+Cj?) 
TP 


VnUn* 


= Re • tr[{3f/3p'l'B-^*3l/3p+^^4>3^/3P*B>VnUn 1 (11) 

where « -(Is-Tt)"^. Note that 
a'B’/ap - 0 
and 

IT . TiPh 

^ - p + i2PH U2) 


where 



Ncx(Nc+M) (Ns+M)x(Nc+M) 



(Ns+M)x(Ns+M) 


matrix A is already available In upper Hessenberg 
form (10) from the computation of the singular 
values. The second Is the computation of the 
eigenvectors and generally involves a low-order 
complex matrix. 


Optimization Schemes 


Let us assume that the minimum over the 
frequency domain of the singular value o^(I+KG) of 
a stable system Is oj^. It is desired to 
Increase to a desired value oq 
illustrated In figure 3. An Increased gjy| 
results In better gain and phase margins of the 
system as shown In equation (7) and figure 2. 
Optimization schemes to achieve this objective 
using the gradient Information of equation (10) 
are described next. 

Unconstrained Minimization Approach 

In the unconstrained minimization approach* a 
single objective function J is minimized by 
changing the design variables p. Since g_(I+K6) 
is less than gj) over a range of frequencies 
Instead of at a single point, all of the 
violations where o^(I+KG)< aj) are represented by 
a single cumulative measure J. 


J(p)-^ (Max{0,[oj)-o(jMi,p)]} )2 (15) 


Equation (11) can be written as 


3on{l+lCG) 

3p 




H»¥vnUn 


(13) 

Using the matrix trace properties Re tr(A)»Re 
tr(A*), a matrix relation for the gradients with 
respect to all of the elements of ? can be written 
as 


3 an(I+KG) ^ __ * 

Re[{Il^(C«t2)*){fljBvnUn*}] (14) 

(Nc+M)x(No+M) 


The summation 1$ taken over a large number of 
frequency points where both the choice of the 
frequency range and spacing of the frequency 
points In a frequency range are left to the 
designer. A geometric description of the 
cumulative objective function Is shown In figure 
3. The objective Is to minimize (preferably 
reduce to zero) the shaded area below the ^ 
line. A conjugate gradient algorithm (11) is used 
to search for the controller design variables p 
which minimize J without allowing g_ to go near 
zero during the search process. Not allowing o to 
go near zero is particularly Important to avol"? 
destabilizing the system during the linear search 
process, especially when g^ has sharp drops at 
specific frequencies. The method Is expected to 
work when g_ and 3o_/3p variations with frequency 
are not too large over small frequency ranges. If 
J can be reduced to zero, then the minimum 
singular value reaches aj) or higher. 


The complex conjugate transpose of equation (14) 
gives equation (10).* The gradient expressions for 
the matrices (I+KG)’^, KG, KG(I+KG)“1 can be 
obtained in the same manner and are given In the 
Appendix. Two additional computations are 
Involved in computing singular value gradients at 
each frequency point. The first computation Is 
the solution of a set of (Ng M) simultaneous 
equations and Is relatively Inexpensive since the 


Constrained Minimization Approach 

In the constrained minimization approach, an 
objective function J Is minimized with respect to 
the design variables p subject to the Inequality 
constraint g<0. In this approach the cumulative 
measure of aTl of the violations g_(I+KG)< oo 
treated as the constraint (12). The objective 
function J and the constraint g are defined as 
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J(p) = 1/2 tr[cTc] 

(16) 

g(P) “X (Max{ 0 , [oj)-a ( ju i, p)] } )2 

(17) 


i 


The choice of J in equation (16) is desirable 
since a lower C is reflected in lower control 
activity. Other choices of J are possible. In 
equation (17), the summation is taken over a large 
number of frequency points as before. A geometric 
description of the cumulative constraint is shown 
in figure 3. The objective is to reduce the 
shaded area to zero by satisfying the inequality 
constraint g£0. Although the present paper is 
confined to a single constraint, additional 
constraints on responses and singular value bounds 
at other points in the loop can be considered for 
an overall design. The method of feasible 
directions (11) is used to search for the 
controller design variables p which minimize J 
subject to g£0. The method uses the objective 
function and constraint gradient information to 
determine a parameter move direction and a scalar 
multiplier in the usable-feasible direction to 
satisfy all constraints. When the constraint 
condition is satisfied, then g<0 which implies 
^>op for all from the definition of g in 
equation (17). 


Numerical Results 


Numerical results are presented for a two- 
input/two-output system which represents a drone 
aircraft with a lateral attitude control system 
(4 and 8). A nominal controller is available for 
comparison. The present method is used to 
increase the robustness by redesigning the nominal 
controller. A block diagram of the drone lateral 
attitude control system is shown in figure 4 
(8). The plant state vector x is defined as 

X = [3 4* !i> 4. 61/20 6 2/20 J 'T 

The plant matrices F, Gy, and H as defined in 
equations (1) and (2) are given In table 1. The 
nominal controller matrices A,B,C, and D as 
defined in equations (3) and (4) are given in 
table 2. The eigenvalues of the nominal open-loop 
and closed-loop system are given in table 3. The 
eigenvalue at X » 0.1889 ± jl.051 results in an 
unstable dutch roll mode. The elements of the 
input vector are the eleven and rudder actuator 
commands, respectively. All gain and phase 
changes are considered at the points X in figure 
4. The minimum singular value of the return 
difference matrix (I+KG) over the operating 
frequency range is plotted in figure 5 for the 
nominal system. The minimum singular value is 
constant at 0.35 over low frequencies, then drops 
to its lowest value of 0.25 near 1.2 radians/sec 
which is close to the frequency of the unstable 
open-loop pole. The minimum singular value 
approaches unity asymptotically as KG attenuates 
at higher frequencies. Using the stability 
condition given in equation (7), the stability is 

guaranteed if cr{L*l-I) < 0.25. This can be 
Interpreted In terms of gain and phase margins 


using figure 2. The guaranteed simultaneous gain 
margins are -2.0 d 8 and 2.5 dB ( 4 > 1 - 4^2 - 0 ). 

The simultaneous phase margins are ±15* (ki ■ 
k 2 - 0 dB). 

Figures 6 a and 6 b show the gradients of 
o_(I+KG) with respect to the nominal controller 
parameters aji, a22****d22* The 
location of these parameters in the block diagram 
is shown in figure 4. The elements bn and 
b 22 do not show up in figure 4 since their unity 
values are embedded in the controller structure. 
The gradients with respect to c^ and dn are 
quite large. The gradients with respect to other 
diagonal elements a^, a? 2 * ^22* ^22* 
etc. are relatively small. These gradients 
attenuate to zero before 10 rad/sec except for the 
one with respect to d]^j which attenuates at 30 
rad/sec. It may be noted that although the off- 
diagonal elements are zero, the singular value 
gradients with respect to them are quite large. 

Results of unconstrained minimization (Design 
1 ) . Unconstrained minimization is performea using 
c\i and d ?2 as the design parameters. The 
desired minimum singular value oj) is 0 . 6 . The 
equality relations dn • 0 and C 22 ■ ^22 • 
d ?2 are maintained to satisfy the 1 /s and 
s/{s+ 2 ) structure of the nominal control law. 

Hence the design parameters are basically 
proportional to gains in each loop. Although the 
convergence pattern of both the objective function 
and the design variables are not shown, the 
objective function reduces to zero in one 
iteration and the values of C 22 3 nd d 22 ^re 
0.13 and 9.69, respectively. Note that C 22 “ 

-2d22 = -19.38. The singular value plot is 
shown in figure 7 as design 1 . The minimum 
singular value om is 0.6 as desired. The price 
paid for the higner value of ojwj is the loss of 
rapid attenuation at higher frequencies. With 
OM » 0 . 6 , the guaranteed gain margins are - 4.1 
dB and 8.0 dB and the phase margins are ±35®. 

This reflects substantial Improvement over the 
nominal stability margins. The eigenvalues of the 
closed-loop system are given in table 3 . 

Results of constrained minimization (Design 

Next the same problem is solved using the 
constrained optimization approach defined as 
design 2 . The desired minimum singular value ^ 
is again 0.6. The objective function is chosen as 
defined by equation (16). The convergence pattern 
for design 2 is shown in figure 8 . The J and g 
are normalized by their starting value Jq and 
gy, respectively. The constraint is satisfied 
in three iterations but at the cost of increased 
J. The values of cn and d 22 after five 
iterations are 2.08x10-^ and 5.91, respectively. 
The corresponding singular value plot is shown in 
figure 7 as design 2 . The minimum singular value 
0 . 68 . The loss of attenuation at higher 
frequencies is much less as compared to design 1 , 
probably because the algorithm tries to minimize 
the growth of 0.S{c\i^ ■** <^ 22 ^) well. 

The eigenvalues of the closed-loop system are 
given in table 3. 
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As a general rule, an Increase In robustness 
Is accompanied by degraded response and Increased 
control activity. This effect Is examined from 
time response plots of the closed-loop system 
using the nominal, design 1, and design 2 
controllers presented In figures 9a through 9e. 

The Input Is a unit ramp-hold eleven conmand which 
rises linearly from 0 to 1 In 0.4 seconds as shown 
In figure 9a. The sideslip 8 response Is shown In 
figure 9b. The Increase In sideslip from the 
nominal are roughly four tiroes for design 1 and 
twice for design 2. Figures 9c and 9d show that 
the roll and yaw rates are 10 to 20% lower than 
nominal. The eleven activity Increases by 25% for 
design 1 and 10% for design 2. The Increase In 
rudder activity Is roughly three times for design 
1 and twice for design 2 with large Initial 
overshoot. 


Conclusions 


A method for designing feedback controllers 
to Increase the robustness of multiloop systems 
has been presented. Gradients of the singular 
values of the return difference matrix with 
respect to design variables In the controller were 
derived analytically. A cumulative measure of the 
singular values and their gradients was used with 
a numerical optimization algorithm to Increase the 
system's robustness. 

A numerical example was given to Illustrate 
the method. For the example, a nominal controller 
was available. The present method was used to 
design a new controller which provided Increased 
robustness. The global minimum singular value was 
Increased substantially using both the 
unconstrained optimization approach and the 
constrained optimization approach. For both of 
these cases, the time response of the system using 
these controllers was degraded. Some high 
frequency attenuation was lost. For a better 
overall design, more constraints need to be added. 
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The singular value gradients of some useful 
atrlces with respect to the controller quadruple 
Is presented In this appendix. The left and 
right eigenvectors Up and Vp In each expres- 
sion belong to the singular value of that 
particular matrix. 


3on(l^KGH ^ .Re[tWa^VnU„*[-(I+KG)-^i (A’D 
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3o„(KG(I+KG)-l . I /« 

- 5 ^: Re[H*a^VnUn*[-(I+KG)"l{ 


®-Tl 


T5»aT]] 


G($) K($) 



where *a = (Is*^ ^ 


(A.4) 


Fig. 1 Block diagram of a multi loop system. 


The gradient expression for a(GK) etc. at the 
output can be derived similarly starting from 
equation 9. 
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Table 1 Plant matrices F, 6y, and H for drone 
lateral attitude control system 



*,08527 

-0.0001423 

-0.9994 

0.04142 

0 

0.1862 


-46.86 

-2.757 

0.3896 

0 

-124.3 

128.6 

F = 

-0.4243 

0 

-0.06224 

1 

-0.06714 

0 

0 

0 

-8.792 

0 

-20.46 

0 


0 

0 

0 

0 

-20.0 

0 


_ 0 

0 

0 

0 

0 

-20.0 


"o o" 







0 0 
0 0 
0 0 
1 0 
0 1 

II 

X 

1 0 

0.07 1 

0 0 
0 0 

:] 



Table 2 Controller quadruple matrices A, B, C, and 
D for drone lateral attitude control 
system 



Fig. 2 Universal diagram for gain-phase margin 
evaluation. 


c.p-“ “ 1 “ 1 

[o -4.U6J LO 2.0Wj 


Table 3 Eigenvalues of drone lateral attitude 
control system 


MODE 

OPEN LOOP 

CLOSED LOOP 

NOMINAL 

DESIGN -1 

DESIGN -2 

1 

-0.03701 

-0.6911 

-0.0386 

-0.01399 

2. 3 

0.18891 ) 1.051 

-0.25531 ) 1.187 

-0.6436 1 ) 0.823 

-0.5336 1 ) 0.959 

4 

-3.25 

-2.60 

-6.225 1 ) 2.342 

-3.866 1 ) 2.276 

5 

-20.0 

-18.70 

-11.11 

-16.06 

6 

-20.0 

-20.15 

-20.02 

-20.00 

7 


0 

0 

0 

8 


-2.261 




SINGULAR 
VALUE 
fl(I + KG) 



FREQUENCY, u) 


Fig. 3 Geometric description of cumulative 
objective function and cumulative 
constraints. 
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Fig. 4 Block diagram of a drone lateral attitude 
control system. 



FREQUENCY cj. rad/sec 

Fig. 6(b) Gradient of singular value a(I+KG) with 
respect to controller parameters cn, 
C22» ^^11* and d22 (nominal). 



FREQUENCY w. rad/sec 

Fig. 5 Maximum and minimum singular value vs. 
frequency for drone attitude control 
system (nominal). 



FREQUENCY w. rad/sec 

Fig 7 Minimum singular value vs. frequency plot 
for nominal and optimized control laws. 



FREQUENCY u. rad/S8C 

Fig. 6(a) Gradient of singular value a(I+KG) with 
respect to controller parameters an, 
®22» ^11 » and b22 (nominal). 



NO. OF ITERATIONS 


Fig. 8 Convergence pattern of constrained 
minimization (Design 2). 
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Fig. 9(a) Ramp- hold eleven cortmand. 


Fig. 9(d) Yaw rate response to ramp-hold 
eleven command. 


r 



Fig. 9(b) Side slip response to ramp-hold elevon 
command. 



Fig. 9(e) Elevon deflection response to ramp-hold 
elevon command. 


ROLL RATE, 
deg/sec 




Fig. 9(c) Roll rate response to ramp-hold 
elevon conmand. 


Fig. 9(f) Rudder deflection response to ramp-hold 
elevon command. 
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